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We hypothesize that in a non-metallic crystalline structure under extreme pressures, atomic wave- 
functions deform to adopt a reduced rotational symmetry consistent with minimizing interstitial 
space in the crystal. We exemplify with a simple numeric variational calculation that yields the 
energy cost of this deformation for Helium to 25%. Balancing this with the free energy gained by 
tighter packing we obtain the pressures required to effect such deformation. The consequent modi¬ 
fication of the structure suggests a decrease in the resistance to tangential stress, and an associated 
decrease of the crystal’s shear modulus. The atomic form factor is also modified. We also compare 
with neutron matter in the interior of compact stars. 


I. INTRODUCTION 


When a solid is compressed, there is great energetic ad¬ 
vantage in reducing the volume occupied per atom. Once 
an optimal close-packing structure is already in place, 
this can be achieved by a reduction of the interatomic 
spacing, or, what we study in this small note, a deforma¬ 
tion of each of the atoms forming the solid. For a simple 
theoretical study Hydrogen is not the most appropriate 
study case because it forms covalent molecules, a four- 
body problem, so we settle to address Helium instead. 
It should be obvious to the reader that the ideas are ex¬ 
tendible to other atomic systems, but the computations 
would be quite more challenging. 

One can think of two obvious modes to deform the 
atomic wavefunction. The first is a simple compressional 
mode where the atom’s radius is reduced but the lattice 
remains the same. The second is a mode in which the 
atom is deformed so as to abandon the spherical symme¬ 
try that it adopts in vacuum and instead makes a transi¬ 
tion to an octahedral (cubic) symmetry consistent with 
the crystal axes, but maintaining the atomic volume con¬ 
stant. Because of the larger novelty of this concept, this 
is the mode that we will be addressing; it is plain that 
both should act on a real system under compression. 
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FIG. 1: Energy cost to deform para-He as function of the 
cubicity p (triangles) and reduction of the lattice volume per 
atom (circles, normalized to the P = 0 Van der Waals radius). 


II. VARIATIONAL DEFORMATION ENERGY 


Deforming the Is wavefunction from spherical symme¬ 
try to cubic symmetry costs energy that we now roughly 
estimate. We employ Jacobi coordinates for the three- 
body nucleus (center of mass) and two electron system, 


P 2 ~Pi 
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7TA = 


Pz-PiT 


7T C 


= PiT 3 - (i) 


Experimentally, Helium crystallizes at T ~ OK under 
a pressure of 2.5 MPa [lj in an hexagonal-compact (hep) 
structure, but presents also a face-centered cubic struc¬ 
ture (fee) with a triple point fluid -fee-hep at T = 15K and 
P = O.IGPa 0. Through X-ray diffraction studies, it has 
been determined that at T = 300 K and high pressures 
between 15.6 and 23.3GPa, He is in the hep phase [3], and 
moreover the crystalline structure is deformed, diminish¬ 
ing the cell parameters due to the decrease in interatomic 
distance (the pure compressional mode). To our knowl¬ 
edge, the second mode of atomic deformation, breaking 
the central symmetry around the nucleus, has not been 
reported in the literature. 


As Rayleigh-Ritz variational ansatz we take 
i>ipi,P2,P3) = with p, A be- 

ing our parameters to minimize the atomic energy for 
each family of wavefunctions. 

We adopt = e~^ Kxx+ ' Kx y +7Txz ) and ijj p = 

_ / T] I 77 T? 

e \ P x^ Py ^ P z) 5 with 7] a “cubicity” parameter. 

For rj = 2 the symmetry of the orbital is spherical. For 
r] > 2 we have a hyperellipsoidal function with octahedral 
symmetry, and for r] -T oc (in practice 77 ~ 10 is very 
large already), the shape is perfectly cubic. 

The binding energy is obtained by numerically solving 
the variational problem for the Hamiltonian of the non- 
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FIG. 2: Pressure needed to deform the Helium atom wave- 
function from spherical to hyperellipsoidal of given cubicity. 


relativistic He atom (neglecting spin and nuclear recoil), 
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with 7^ = 5. Laboratory studies have reached pressures up 
to 80 GPa with He and the full 500 GPa with Hidrogen |T] 
with the diamond anvil method, so that the wavefunction 
deformation that we propose, while as yet unseen, is not 
very far from sight. 

Atomic orbital deformation is a concept related to lo¬ 
calized electrons around a nucleus. In a metal, delocal¬ 
ized electrons respond to very different physics. So if 
metallization would occur at a pressure less than what 
we estimate for cubicity, there would be no hope. 

For He however there are sophisticated density- 
functional theory computations in the Generalized Gra¬ 
dient Approximation (GGA) [5!, that suggest that the 
gap to the conduction band at T = OK vanishes at 
p = 17.4 g/cm 3 , equivalent to a pressure of P= 17.0 TPa. 
Alternative simulations based on Diffusion Monte Carlo 
(DMC)|6], yield an even larger density of 21.3(1) g/cm 3 , 
corresponding to a pressure P = 25.7 TPa. Both these 
figures are well above the 0(l)TPa at which a local, cu¬ 
bic wavefunction in the insulating phase would occur as 
suggested by our estimates. 


IV. DECREASED SHEAR MODULUS 


H 7 / = Ei/j 


mm 

(V#) 


The 6 and 9-dimensional integrals resulting are computed 
with Vegas [4], a standard Monte Carlo method. The 
computed energy for this interpolating family of func¬ 
tions is shown in figure [l] It costs about 9-10 eV to to¬ 
tally deform the Helium atom to cubic symmetry; the 
accuracy of this computation is about 25% from the 
simple-minded variational approach and ignoring the in¬ 
teratomic inter acti on (appearing only through the pres¬ 
sure in section [TTT| ) . For 77 = 2 the binding of He is un¬ 
derestimated to be 59 eV instead of 79 eV, but since we 
subtract E(rj) — E( 2), a good part of the error will cancel 
out in the difference. We feel that higher accuracy at this 
stage, without experimental data, would be meaningless. 


III. COMPUTATION OF THE PRESSURE 
REQUIRED FOR A GIVEN DEFORMATION 

The deformation energy can be provided by 
Helmholtz’s free energy A E = P AV, obtained fol¬ 
lowing the reduction of the total volume of the hep 
lattice structure since interstices are reduced. The 
computation of AV is a geometric problem that we 
defer to the appendix. The outcome has been plotted in 
figure [T| as function of the cubicity. 

The pressure for which the Helmholtz free energy bal¬ 
ances the wavefunction deformation cost is plot, as func¬ 
tion of the interpolating cubicity parameter, in figure |2j 
In the limit of total deformation 77 —» 00 , P = 1.17 TPa; 
for P = 500 GPa the deformation is already significant, 
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FIG. 3: Top: Representation of the atomic displacement. 
Bottom: Variation of the displacement angle with cubicity 
when the crystal undergoes shear stress. 


When optimally packing spheres, the successive crys- 
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tal planes intertwine. A lateral displacement requires 
(classically) lifting a plane so the spheres can slide (we 
neglect tunnelling, possible since the atoms are not per¬ 
fectly rigid bodies, but that opposes Pauli’s principle). 
This lifting requires a force to compensate the pressure. 
But at large P, the atomic deformation allows optimal 
packing without atoms intruding in neighboring planes 
(controlled by the function g(rj) in Eq. (A3)), since the 
necessary 2 R(rj) separation between plane centers is con¬ 
sistent with optimal packing for large g. Thus, we ex¬ 
pect the shear modulus to decrease. A calculation of the 
complete shear modulus requires also the weak He inter¬ 
atomic potential; we will refrain from addressing it and 
only consider differences in the shear modulus between 
deformed and undeformed He. 
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top panel). This is in turn a function of g through the 
functions P, / and g in the appendix. Numerically it 
falls quickly with g (figure [3j bottom panel). When it 
reaches zero the planes can slide sideways without resis¬ 
tance from the pressure, only cohesive forces retain them 
(taking them into account adds a term to Eq. ^ and our 
actual approximation is to take this independent of g). 

The force exerted by the pressure that resists plane 
lifting 1 * is obtained from energy balance, F L = PAH, 
while the volume increase is 

AV = V{L)-V(L = 0) = A(L+d)-Ad = AL (4) 

(d is the plane separation in the undeformed structure). 
Thus, the contribution of the geometry distortion to G is 

AG = G V - G 2 = P [tan (/^(L)) - tan (/3 2 (L))] (5) 

The outcome is plotted in figure [4j top plot. As the 
cubicity g increases, the shear modulus G decreases, so 
| AG| increases to a saturation of about 365 GPa, or about 
25% of the pressure needed to attain that deformation. 
The actual observable prediction is that the shear mod¬ 
ulus flattens out at large pressure, instead of growing 
linearly with it. 

Finally, we provide a check that we remain in the linear 
elasticity regime in figure[4j bottom panel, where we show 
that the dependence of G with the vertical displacement 
L is very small and only appreciable for small cubicities, 
and even then it is a correction less than 2 GPa. 


V. FORM FACTOR 

A direct way of ascertaining the structure of a charge 
distribution is to probe it by scattering an electron or 
photon beam. Shape information about the target is en¬ 
coded in the elastic form factor. This is defined so that 
the charge density of the object under study is normal¬ 
ized to 1 , f d 3 xp(x ) = 1. The charge density is easily 
computed from the variational wavefunctions ijj\ and Vv 
as in section |II] (do not confuse this charge density with 
Jacobi’s coordinate for the three body problem). Since 
we aim to calculate the form factor of an atom, the non¬ 
relativist ic approximation is sufficient, so that 


FIG. 4: Top: Absolute value of the decrease of the shear 
modulus. Bottom: Near independence of AG with the plane¬ 
lifting L for two cubicities, as a check of linear elasticity. 

The shear modulus G is the ratio between stress r = 
F/A and deformation e = x/L. The displacement paral¬ 
lel to the crystal planes is x, and L their lifting, so 

5 A x/L l/tan/3 

as expressed in terms of /?, the minimum angle of at¬ 
tack that allows tangential displacement (see figure [3j 


F (r ! ) = / (|y(6) 

with F( 0) = 1 due to the charge normalization of the two- 
electron cloud (overall neutrality is guaranteed by the He 
nucleus, which, being approximately point like, provides a 


1 It is of notice that a shear force produces a small crystal dilata¬ 

tion, due to the normal component of the force in the decompo¬ 

sition at the tangent point between atoms, in a similar way to 
how a body ascends an inclined plane when pulled laterally. 
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FIG. 5: Form factor for different cubicities (grade of the hy- 
perellipsoidal wavefunction) as function of the modulus of the 
incident photon’s momentum. 

trivial constant contribution that we understand as sub¬ 
tracted). 

We carry out this integration numerically. We imagine 
that the compressed atomic system (whose deformation 
we should like to measure) is analyzed with X-ray pho¬ 
tons of momentum q. Such X-ray diffraction experiments 
have been carried out as already discussed [3], though at 
smaller pressures as we would need here, but they might 
be performed in the future. 

We now briefly observe what simple changes of the 
form factor one would hope to find due to the atomic 
deformation. 

The numeric computation of the form factor as func¬ 
tion of ||g|| is given in figure |5j Since our wavefunction 
ansatz is a simple exponential, the Fourier transform 
yields a falling rational form factor as seen in the fig¬ 
ure. The difference in the momentum dependence to the 
spherical case is subtle and likely not isolatable from ex¬ 
perimental data with finite error bars. We observe that 
the form factor decreases a bit quicker for higher defor¬ 
mation. 

More interesting is the breakup of azimuthal symmetry 
by the cubic deformation. We give the change of the form 
factor with that atomic deformation as function of the po¬ 
lar scattering angle in figure [ 6 ] (top panel). The form fac¬ 
tor changes significantly with varying cubicity; for g = 2 
the form factor is independent of the measurement angle 
because the charge distribution is spherically symmet¬ 
ric, property maintained by the Fourier transform (top 
squares following a horizontal line) while for deformed 
wavefunctions the form factor is not only smaller but 
also acquires a remarkable angular dependence which is 
the tell-tale signature. 

There is of course also a dependence on the azimuthal 
angle that we show in figure [ 6 ] (bottom panel) for a \\q\\ 
fixed and varying polar angle. The reduced rotational 
symmetry of the deformed atoms for a fixed cubicity pa¬ 
rameter causes some degeneracies that provide a practical 
check of the computer code. 
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FIG. 6: Top: Form factor as function of 0 q for selected cubic¬ 
ities, with fixed ||g||=4eV and cj) q — 'K/2. Bottom: Variation 
with the polar angle 6 q of the form factor for fixed ||g|| =2 eV 
and rj = 6 showing the dependence on the azimuthal angle. 
Note that pairs of numerical data sets fall on top of each other 
due to the 7r/2 rotational symmetry of the hyperellipsoidal 
wavefunctions. 


VI. OUTREACH AND DISCUSSION 

In section [IT] we employed a direct variational calcu¬ 
lation, but there is another way of thinking of a wave- 
function with cubic symmetry, in terms of a spherical- 
harmonic decomposition: a linear combination of the Is 
orbital (with weight ~ 74%), a 3d orbital (about 25%) 
and residual amounts of higher g , i, etc. orbitals. The 
energy estimated this way is similar to what we have ob¬ 
tained in section El 

Extreme pressures are reached not only in the labora¬ 
tory but also in compact stars. In the cooling of white 
dwarves, the metal-insulator transition of He changes the 
carriers of energy from e~ to 7 0. Since metallization 
occurs deep in the star B Et there is ample room for 
atomic deformation in the crust, and as we suggest it 
affects the mechanical properties of the star. 

Neutrons in neutron stars are also subject to extreme 
pressures, and though often treated as a superfluid, they 
may crystallize uni] Heisenberg’s uncertainty principle 
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FIG. 7: Absolute value of the decrease of the shear modulus as 
a consequence of neutron deformation in a neutron star, were 
an anisotropic crystal to develop, as function of the pressure. 
The cubicity is in the interval (2,18). 


AxAp ~ 1 A E = ~ 2 m^Ax 2 su gg ests the Pres¬ 

sure necessary to confine the neutron in a lattice cell, with 
the volume gain estimated from the star’s density and 
Kepler’s packing fraction M ■ Ar— (/crystal /fluid) Veil 
with Veil rsj 2i7 3 , R n = 0.8775 fm. In appropriate units 
( m k = 138 MeV), the pressure is P ^ 1.5m^, which oc¬ 
curs uni at a depth of 3 km below the star surface (typical 
radii are 10 km). 

In such conditions, neutron deformation is a possibility, 
and we calculated in mi that the neutron mass as a 3- 
quark system increases by 150 MeV between 77 = 2 and 
77 = 18. The necessary pressure saturates at about ^ 
llm^, and the change of the shear modulus as function 
of the pressure is shown in figure [7] this may be of interest 
to compute the sismic response of neeutron stars m- 

To conclude, we feel that wavefunction deformation is 
an interesting qualitative phenomenon possible in highly 
compressed systems, with applications beyond labora¬ 
tory experiments to astrophysics, and hope that it is born 
in future, more sophisticated calculations, as well as in 
experimental work. 


Appendix A: Packing of hyperellipsoids 


FIG. 8: To keep the hyperellipsoid’s volume constant notwith¬ 
standing the deformation, its (short) radius decreases as 
shown. 



FIG. 9: Atomic disposition in a unit cell at close-packing. 


for the volume m , 


V{rj) = 8R(r]) 3 


r( 1 + - 

77 


r 



(Al) 


If we understand the radius of the hyperellipsoid as the 
minimum distance from its surface to its symmetry cen¬ 
ter, the relation between this radius an d th at of the 
sphere of equal volume is, because of Eq. (Al), 


Hyperellipsoidal coordinates, </>E [0, 27 t], #E [0 ,7 t], with 

x = ai cos 2 / 77 # cos 2 / 77 ^, y = <22 cos 2 / 77 # sin 2 / 77 ^, and z = 

as sin 2 / 77 # are used for the hyperellipsoid’s volume. 

First we obtain the area of the hyperellipse obtained 

upon cutting the figure with a horizontal plane by means 

of Green’s theorem A = \ <f> c (xdy — y dx ), which is 

A = — b(~, - + l] in terms of Euler’s beta function. 

7] y T] 5 7] J 

Letting now the axes of the hyperellipse vary with height 
or polar angle, a(#), #(#), we obtain the known formula 


R(v) = 


"7T1 V 3 
. 6 . 


i + 5 

T/ 


1 1/3 


i + l 

77 


Ro 


(A2) 


The resulting R(r]) is pictured in figure (18]). Obviously 
for large cubicity, lim^oo 77 ( 77 ) = [7r/6p 3 i7o because 
the cube’s volume is 8R(rj oo) 3 . 

We are ready to compute the packing fraction; consider 
the first octant of a given hyperellipsoid I and the adja- 
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cent ones in a close-packing structure that are tangent in 
that octant, numbered II through IV in figure [9| 

When 77 = 2, the sphere’s centers are (0, 0, 0), 
(2Ro, 0 , 0 ), (i?o? VSRo, 0 ) and (Rq, y/l/3Ro, ^/S/3 Rq). 
Upon deforming the figures, they may move closer 
since space is used more efficiently, and the centers 
become ( 0 , 0 , 0 ), (2R(rj), 0 , 0 ), ( R(rj ), f(rj)R(r 7 ), 0 ) and 
(x, 7 /, g(rj)R(r 7 )), where the auxiliary functions R(rj), (de¬ 
termined numerically), x, y , f(rj) = (2 V — 1 ) 1//r? and 77 ( 77 ) 




2 "-l- 


/ / Q ) 2 -1 

V 2/(7?) 


? 7 -] 1/77 


(A3) 


are found requesting that the four hyperellipsoids 


x V +y V + z V = R( v )r> (A4) 
(x - 277(7?))" + y" + = 77 ( 7 ?)" (A5) 

(x - 7?(t?))" + (t? - /(t?)7?(t?))" + 2 " = 7?(t?)" (A 6 ) 

(^-7?(7?))"+(^^77(7?))% 

(z - 3 (t?)7?(t?))" = 77(7?)" (A7) 


are tangent at precisely one point. 

Once we have the size and position of the deformed 
figures we can calculate the occupied volume fraction. 
Without loss of generality we may choose Rq = 1 and take 
as total volume a cube of radius R{rj) with a vertex on 
the origin, that just leaves hyperellipsoid II out; but III 
and IV, because of the closely packed structure, fill some 
of the interstitial space. The octant of hyperellipsoid I 
occupies Vi «■? ^. The other two occupy 


Vru = 


V w = 
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/' 

Jfl 


R-^-lR-fRl 71 ) 1 2 3 4 5 6 7 ^ JfR-{R^-\x-R\ T1 ) 1 /' n 
( R ri-\x-R\«-\y-fR\n)Vrt 

/ dxdydz (A 8 ) 

Jo 


L 

l 


f 

1,r > Jy 

gR-(^RV-\ X -R\V-^y-£^zl R ^ 


iR-iRv-lR-gRlvy/n Jy inf 
r R 


(A9) 


^ i dxdydz . 


The resulting occupied fraction can be subtracted from 
1 to obtain the wasted or interstitial fraction depicted 
in figure [To| For 77 = 2 we reproduce Kepler’s result for 
optimal packing, about 26% inefficient JJ. As 77 grows we 
see that this number decreases and eventually vanishes 
(the packing of cubes leaves no interstitial space). For 
comparison we also show the rather trivial result for a 
simple-cubic lattice. For qualitative reasoning, we see 
that the actual type of lattice is not very important. 
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FIG. 10: Interstitial volume in the hep and cs structures. 
The volume gained upon occupation of the interstices does 
not differ so much between the two packings. 
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